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Kolmogorov spectrum of superfluid turbulence: numerical analysis of the 
Gross-Pitaevskii equation with the small scale dissipation 
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The energy spectrum of superfluid turbulence is studied numerically by solving the Gross- 
Pitaevskii equation. We introduce the dissipation term which works only in the scale smaller than 
the healing length, to remove short wavelength excitations which may hinder the cascade process 
of quantized vortices in the inertial range. The obtained energy spectrum is consistent with the 
Kolmogorov law. 
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The physics of quantized vortices in liquid 4 He is one of 
the most important topics in low temperature physics [l| . 
Liquid 4 He enters the superfluid state at 2.17 K. Below 
this temperature, the hydrodynamics is usually described 
using the two-fluid model in which the system consists 
of inviscid superfluid and viscous normal fluid. Early 
experimental works on the subject focused on thermal 
counterflow in which the normal fluid flowed in the op- 
posite direction to the superfluid flow. This flow is driven 
by the injected heat current, and it was found that the 
superflow becomes dissipative when the relative velocity 
between two fluids exceeds a critical value p| . Feynman 
proposed that this is a superfluid turbulent state consist- 
ing of a tangle of quantized vortices [2| . Vinen confirmed 
Feynman's picture experimentally by showing that the 
dissipation comes from the mutual friction between vor- 
tices and the normal flow |3(- After that, many studies on 
superfluid turbulence (ST) have been devoted to thermal 
counterflow However, as thermal counterflow has no 
analogy with conventional fluid dynamics, we have not 
understood the relation between ST and classical turbu- 
lence (CT) . With such a background, experiments of ST 
that do not include thermal counterflow have recently ap- 
peared. Such studies represent a new stage in studies of 
ST. Maurer et al. studied a ST that was driven by two 
counter-rotating disks [f| in superfluid 4 He at temper- 
atures above 1.4 K, and they obtained the Kolmogorov 
law. Stalp et al. observed the decay of grid turbulence [|J 
in superfluid He at temperatures above 1 K, and their 
results were consistent with the Kolmogorov law. The 
Kolmogorov law is one of the most important statistical 
laws of fully developed CT, so these experiments show 
a similarity between ST and CT. This can be understood 
using the idea that the superfluid and the normal fluid 
are likely to be coupled together by the mutual friction 
between them and thus to behave like a conventional fluid 
@- Since the normal fluid is negligible at very low tem- 
peratures, an important question arises: even without 
the normal fluid, is ST still similar to CT or not? 

To address this question, we consider the statistical 
law of CT 0. The steady state for fully developed tur- 
bulence of an incompressible classical fluid follows the 
Kolmogorov law for the energy spectrum. The energy is 



injected into the fluid at some large scales in the energy- 
containing range. This energy is transferred in the iner- 
tial range from large to small scales without been dissi- 
pated. The inertial range is believed to be sustained by 
the self-similar Richardson cascade in which large eddies 
are broken up to smaller ones, having the Kolmogorov 
law 



E(k) = Ce 2 / 3 k~ 5 / 3 . 



(1) 



Here the energy spectrum E(k) is defined as E = 
J dk E(k), where E is the kinetic energy per unit mass 
and k is the wave number from the Fourier transfor- 
mation of the velocity field. The energy transferred to 
smaller scales in the energy-dissipative range is dissipated 
by the viscosity with the dissipation rate, which is iden- 
tical with the energy flux e of Eq. |QJ in the inertial 
range. The Kolmogorov constant C is a dimensionless 
parameter of order unity. 

In CT, the Richardson cascade is not completely un- 
derstood, because it is impossible to definitely identify 
each eddy. In contrast, quantized vortices in super- 
fluid are definite and stable topological defects. A Bose- 
Einstein condensed system yields a macroscopic wave 
function &(x,t) = yf p(x, i)e 1 ^ x ' t ', whose dynamics is 
governed by the Gross-Pitaevskii (GP) equation 0, 
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where /i is the chemical potential, g is the coupling con- 
stant, p(x, t) is the condensate density, and v(x,t) is 
the superfluid velocity given by v(x,t) — 2V(j)(x,t), 
The vorticity rot v{x, t) vanishes everywhere in a single- 
connected region of the fluid; any rotational flow is car- 
ried only by a quantized vortex in the core of which 
<fr(x,t) vanishes so that the circulation is quantized by 
4tt. The vortex core size is given by the healing length 
£ = 1/ >f9P- Although quantized vortices at finite tem- 
peratures can decay through mutual friction with the 
normal fluid, at very low temperatures the vortices can 
decay only by two mechanisms; one is the sound emis- 
sion through vortex reconnections [ll|, and the other is 
that vortices reduce through the Richardson cascade pro- 
cess and eventually change to elementary excitations at 
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the healing length scale. In any case, dissipation occurs 
only at scales below the healing length. Therefore, ST 
at very low temperatures, consisting of such quantized 
vortices, can be a prototype to study the inertial range, 
the Kolmogorov law and the Richardson cascade. 

There are two kinds of formulation to study the dy- 
namics of quantized vortices; one is the vortex-filament 
model 12], and the other the GP model. By using the 
vortex-filament model with no normal fluid component, 
Araki et al. studied numerically a vortex tangle starting 
from a Taylor-Green flow, thus obtaining the energy spec- 
trum consistent with the Kolmogorov law |13| . By elim- 
inating the smallest vortices whose size is comparable to 
the numerical space resolution, they introduced the dis- 
sipation into the system. Nore et al. used the GP equa- 
tion to numerically study the energy spectrum of ST |14| . 
The kinetic energy consists of a compressible part due 
to sound waves and an incompressible part coming from 
quantized vortices. Excitations of wavelength less than 
the healing length are created through vortex reconnec- 
tions or through the disappearance of small vortex loops 
[Til fl5| , so that the incompressible kinetic energy trans- 
forms into compressible kinetic energy while conserving 
the total energy. The spectrum of the incompressible 
kinetic energy is temporarily consistent with the Kol- 
mogorov law. However, the consistency becomes weak in 
the late stage when many sound waves created through 
those processes hinder the cascade process 0, 0, . 

Our approach here is to introduce a dissipation term 
that works only on scales smaller than the healing length 
£. This dissipation removes not vortices but short wave- 
length excitations, thus preventing the excitation energy 
from transforming back to vortices. Compared to the 
usual GP model, this approach enables us to more clearly 
study the Kolmogorov law. 

To solve the GP equation numerically with high accu- 
racy, we use the Fourier spectral method in space with 
periodic boundary condition in a box with spatial reso- 
lution containing 256 3 grid points. We solve the Fourier 
transformed GP equation 

ijU(fc,t) = ]fc 2 -M]*(M) 



+^ E £(*i,*)£*(*2,t) 

x$(fc -ki + k 2 ,t), 



(3) 



where V is the volume of the system and $(fc,t) 
is the spatial Fourier component of &(x,t) with the 
wave number fc, numerically given by a Fast-Fourier- 
Transformation [l^. We consider the case of g = 1. 
For numerical parameters, we used a spatial resolution 
Ax = 0.125 and V = 32 3 , where the length scale is 
normalized by the healing length £. With this choice, 
Afc = 27r/32. Numerical time evolution was given by 
the Runge-Kutta-Verner method with the time resolu- 
tion A* = i x icr 4 . 

To obtain a turbulent state, we start from an initial 



configuration in which the condensate density po is uni- 
form and the phase 4>o(x) has a random spatial distri- 
bution. Here the random phase (f>o(x) is made by plac- 
ing random numbers between — 7r to it at every distance 
A = 4 and connecting them smoothly (Fig. The ini- 
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FIG. 1: How to make the random phase <f>o(x) with A = 4 (a) 
and one example of an initial random phase in a cross section 
(/>o{x,y,0) (b). 



tial velocity v(x,t = 0) = 2\7<fio(x) given by the initial 
random phase is random, hence the initial wave func- 
tion is dynamically unstable and soon produces homo- 
geneous and isotropic turbulence with many quantized 
vortex loops. 

To confirm the accuracy of our simulation, we calcu- 
late the total energy E, the interaction energy -ffint, the 
quantum energy E q and the kinetic energy -Ekin |l4| : 



E = 

Eint = 

E q = 

Ev\„ — 



fdx<S>*{-\7 2 +.g/2|$| 2 ]$ 
J dx p 

g jdx\^ 

2 Jdxp ' 
Jdx [V|$|]2 

Jdxp 
jdx [|$|V0] 2 
J dx p 



(4) 



and compare them with those given by the dynamics with 
the different time resolution At — 2 x 10~ 5 or different 
spatial resolutions 512 3 grids. For the comparison, we es- 
timate the relative errors Fi 2 {A) = \ ((A) 1 - (A) 2 )/(A) 1 \ 
and F 13 (A) - \{(A) 1 - (A) 3 )/(A) 1 \ at t = 12 for A = E, 
E- mt , E q and £ki n , where () 2 and () 3 are respec- 
tively the values given by three simulations with differ- 
ent resolutions (i) At = 1 x 1CP 4 and 256 3 grids, (ii) 
At = 2 x 1CT 5 and 256 3 grids and (iii) At = 1 x 1CT 4 
and 512 3 grids. Comparison is shown in TableHJ the rela- 
tive errors are extremely small, which allows us to use the 
resolution (i). Furthermore, the total energy is conserved 
by the accuracy of 1CP 10 in our simulation. 

We introduced a dissipation term in Eq. © to remove 
excitations of wavelength shorter than the healing length. 
The imaginary number i in the left-hand side of Eq. @ 
was replaced by [i — j(k)], where -f(k) = jo6(k — 27r/£) 
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TABLE I: 

-Bkin- 



Dependence of J12 and F13 on E, -Eint, E q and 



with 9(x) being the step function. The effect of this 
dissipation is shown in Fig. [3 We divide up to the 
kinetic energy i^kin into the compressible part E^ in = 
J dx [(|$|V<^>) C ] 2 / / da; p, due to sound waves, and the 
incompressible part E^ — J dx (f)) 1 ] 2 / J dx p, due 

to vortices, where rot(|$|V^)) c = and div(|$|V0)' = 
0,0. Figure |21 shows the time development of E, E^ in , 
E^ in and E^ in in the case of 70 = (a) and 70 = 1 (b). 
Without dissipation, the compressible kinetic energy E^ in 
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FIG. 2: Time development of E, E^n, E^. in and E kin . (a) 
Case with 70 = (b) case with 70 = 1. 



is increased in spite of conserving the total energy (Fig. 
ETa )), w hich is consistent with the simulation by Nore et 
al. [l4j . The dissipation suppresses E£ in and thus causes 
E^ in to be dominant. This dissipation term works as vis- 
cosity only at the small scale. At other scales, which 
are equivalent to the inertial range Ak < k < 27r/£, the 
energy does not dissipate. 

We calculated the spectrum of the incompressible ki- 
netic energy E^ in (k) defined as E^ in — J dk E^ in (k). Ini- 
tially, the spectrum i£ kln (fc) significantly deviates from 
the Kolmogorov power-law, however, the spectrum ap- 
proaches a power-law as the turbulence develops. We 
assumed that the spectrum E^ in (k) is proportional to 
k~ v in the inertial range Ak < k < 2tt/£, and we de- 
termined the exponent 77. The time development of 77 
is shown in Fig. [3] (a). This figure shows that our ST 



satisfies the Kolmogorov law E^ in (k) oc fc~ 5 / 3 for times 
4 < t < 10, when the system may be almost homogeneous 
and isotropic turbulence. We also calculated the energy 
dissipation rate e = —dE^/dt and compared the results 
quantitatively with the Kolmogorov law. Figure [3] (b) 
shows that e is almost constant in the period 4 < t < 10, 
which means that the Kolmogorov spectrum e 2 / 3 k~ v is 
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FIG. 3: Time dependence of the exponent of the spectrum 
77 (a) and the dissipation rate e. (a) 77 for 70 = 1. (b) e for 
70 = 1. (c) r\ for 70 = 0. (d) e for 70 = 0. In (a) and (c), the 
error bars are the standard deviation of the data from the fit. 
The line r\ — 5/3 is shown to compare the results with the 
Kolmogorov law. 



also constant in the period. Without the dissipation term 
j(k), the ST could not satisfy the Kolmogorov law. This 
means that dissipating short wavelength excitations are 
essential for satisfying the Kolmogorov law. Figure [3(c) 
shows the time development of the exponent 77 for 70 = 0. 
In this case, 77 is consistent with the Kolmogorov law only 
in the short period 4 < t < 7; active sound waves af- 
fect the dynamics of quantized vortices for t > 7. Such 
an effect of sound waves is also clearly shown in Fig. 
01 (d). Compared with Fig. 03 (b), e is unsteady and 
smaller than that of 70 = 1 and even becomes negative 
at times. When e is negative, the energy flows backward 
from sound waves to vortices, which may prevent the en- 
ergy spectrum from satisfying the Kolmogorov law for 
t > 7. According to Fig. 01(a), the agreement between 
the energy spectrum and the Kolmogorov law becomes 
weak at a later stage t > 10, which may be attributable 
to the following reasons. In the period 4 < t < 10, the 
energy spectrum agrees with the Kolmogorov law, which 
may support that the Richardson cascade process works 
in the system. The dissipation is caused mainly by re- 
moving short wavelength excitations emitted at vortex 
reconnections. However, the system at the late stage 
t > 10 has only small vortices after the Richardson cas- 
cade process, being no longer turbulent. The energy spec- 
trum, therefore, disagrees with the Kolmogorov law of 
77 = 5/3. And emissions of excitations through vortex 
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reconnections hardly occurs, which reduces greatly the 
energy dissipation rate e. 

Spatially, the system appears to develop full turbu- 
lence. Shown in Fig. 0] (a) is the spatial distribution 
of vortices. The vorticity plot suggests fully developed 
turbulence. The energy spectrum in Fig. 01 (b) agrees 
quantitatively with the Kolmogorov law. We thus con- 
clude that the incompressible kinetic energy of ST with- 
out the effect of short wavelength excitations satisfies the 
Kolmogorov law. This inertial range 0.20 < k < 6.3 is 
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FIG. 4: Vorticity rotv(x,t) and the energy spectrum at 
t = 6.0. The contour in (a) is 98% of maximum vorticity. 
The energy spectrum was obtained by making an ensemble 
average for 20 different initial states The solid line refers to 
the Kolmogorov law in Eq. (0. 



larger than those in the simulations by Araki et al. 01 
and Nore et al. ^ij. The Kolmogorov constant is es- 
timated to be C ~ I0~ - 5 ~ 0.32. This Kolmogorov 
constant is small because this turbulence is dissipative, 
being expected to become larger for a steady turbulence 
with injection of energy. 



In conclusion, we investigate the energy spectrum of 
ST by the numerical simulation of the GP equation. Sup- 
pressing the compressible sound waves which are created 
through vortex reconnections, we can make the quantita- 
tive agreement between the spectrum of incompressible 
kinetic energy and Kolmogorov law. Our next subjects 
are making the steady turbulence by introducing injec- 
tion of the energy and broadening the inertial range of 
the energy spectrum. We will report on this results in 
near future. 
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